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Planning is underway for several space-borne gravitational wave observatories to be built in the 
next ten to twenty years. Realistic and efficient forward modeling will play a key role in the de- 
sign and operation of these observatories. Space-borne interferometric gravitational wave detectors 
operate very differently from their ground based counterparts. Complex orbital motion, virtual inter- 
ferometry, and finite size effects complicate the description of space-based systems, while nonlinear 
control systems complicate the description of ground based systems. Here we explore the forward 
modeling of space-based gravitational wave detectors and introduce an adiabatic approximation to 
the detector response that significantly extends the range of the standard low frequency approxi- 
mation. The adiabatic approximation will aid in the development of data analysis techniques, and 
improve the modeling of astrophysical parameter extraction. 



I. INTRODUCTION 

Gravitational wave astronomy can be broadly divided 
into high and low frequency bands, with the dividing 
line near one Hertz. Seismic and gravity gradient noise 
prevent ground based detectors from exploring the low 
frequency portion of the spectrum, making this source- 
rich region the sole preserve of space-based observatories. 

Ground and space-based interferometric gravitational 
wave detectors operate according to the same general 
principles, but differ in their implementation. Ground- 
based detectors, such as the Laser Interferometer Grav- 
itational Wave Observatory (LIGO) 0, operate in the 
low frequency limit, where the wavelength of the gravita- 
tional waves is considerably larger than the size of the de- 
tector, and most sources are only in-band for a fraction of 
a second. These considerations simplify the description of 
the detector response, which may be well approximated 
by a quadrupole antenna moving at constant velocity 
with respect to the gravitational wave source. However, 
ground-based interferometers employ quasi-fixed rather 
than freely moving test masses, and the output of the 
detector is given by the response of the control loop used 
to keep the interferometer on a dark fringe. This com- 
plicates forward modeling efforts for ground-based detec- 
tors as it makes the detector response non-linear. The 
situation with space-borne detectors is completely the 
opposite. Space-based detectors, such as thcproposed 
Laser Interferometer Space Antenna (LISA) Q, will be 
able to detect gravitational waves with wavelengths that 
range from many times larger than the interferometer to 
many times smaller, and most sources will be in-band 
for months or years, so that the detector's orbital mo- 
tion will impart amplitude, frequency, and phase mod- 
ulations. These effects give rise to a complicated, time 
dependent detector response function Space-borne 
detectors typically have large arm-lengths (5 x 10^ m for 
LISA) that vary with time, which prevents them from 
operating as traditional interferometers. Instead, the in- 
terferometer signals are produce in software from phase 
differences measured in the detector using a procedure 
know as Time Delay Interferometry (TDI) |^. Despite 



these complications, the detector response remains lin- 
ear, which greatly simplifies forward modeling efforts. 

Forward modeling plays a key role in the design of 
any new scientific instrument, and is especially impor- 
tant when the instrument is the first of its kind. Work 
is now underway toproduce an end-to-end model of the 
LISA observatory 0. Key ingredients include accurate 
modeling of the spacecraft orbits and photon trajectories 
(this includes the effects of gravitational waves), realistic 
simulations of the time delay interferometry used to can- 
cel laser phase noise, and experimental characterization 
of the various noise contributions. A good end-to-end 
model can help to make design trade-offs, and to avoid 
costly mistakes. Forward modeling can also be used to 
develop and test data analysis strategies. While we focus 
our attention on LISA, our forward model can be used to 
study other proposals for space-borne gravitational wave 
detectors, such as the Big Bang Observatory 0. 

Work on various elements of the LISA end-to-end 
model have been under development for some time. Mod- 
eling of the detector response has its roots in the Doppler 
tracking of spacecraft Q. Results were initially derived 
for a static array with equal arm- lengths |^ . Follow- 
ing the discovery of Time Delay Interferometry |^ , these 
results were extended to a static array with unequal arm- 
lengths 1^, 0, ^2 ■ The orbital motion of the array was 
first incorporated in the low frequency limit Il3|. and 
later extended to the full detector response With 
the full response function in hand, we have developed an 
open source software package called The LISA Simula- 
tor that takes as its input an arbitrary gravitational 
wave and returns as its output the simulated response of 
the LISA observatory. The main purpose of The LISA 
Simulator is to aid in the development of data analysis 
tools 0, 0, [Til ITtI , but its modular design allows it to 
be extended into a full end-to-end model. For example, 
the static modeling jl^J of the TDI implementation could 
be incorporated into The LISA Simulator, as could more 
realistic spacecraft orbits and experimentally determined 
noise spectra. 

The value of a realistic end-to-end model has already 
become apparent with the discovery of flaws in the ini- 
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tial TDI scheme caused by the rotation of the array [T^ , 
time dependence of the arm-lengths [23|, and prob- 
lems with clock synchronization in a moving array |2ll |. 
These difficulties require modification of the TDI vari- 
ables Ulll^,!!! and/or changes in the mission design. 

On the other hand, a highly realistic end-to-end sim- 
ulation necessarily consumes a great deal of computer 
resources, and delivers a fidelity that exceeds the re- 
quirements of many data analysis efforts. Indeed, when 
searching a large parameter space, fidelity must be sac- 
rificed in favor of speed. To this end we have developed 
an approximation to the full LISA response that extends 
the low frequency approximation by two decades. The 
motion of the array is stroboscopically rendered into a 
sequence of stationary states, yielding an adiabatic ap- 
proximation to the full response. The adiabatic approxi- 
mation allows us to write down a simple analytic expres- 
sion for the response function in a mixed time/frequency 
representation. For sources with a few dominant harmon- 
ics, such as low eccentricity, low spin binary systems at 
second post-Newtonian oder, the adiabatic approxima- 
tion provides a fast and accurate method for calculating 
the LISA response. 

The outline of this paper is as follows: In Sec^Jwe de- 
scribe the orbits of the interferometer constellation and 
describe how various effects enter into the detector re- 
sponse. In Sec. lIIII we review the expression for the com- 
plete response of a space-borne detector. (An alternative 
derivation of the full response is given in Appendix B). 
In Sec. llVl we show some applications of the general for- 
malism using The LISA Simulator. In Sec. we explore 
the limitations of the low frequency approximation, and 
in Sec. I VII we introduce the adiabatic approximation and 
demonstrate its utility. We finish with an application, us- 
ing the adiabatic approximation to determine when LISA 
can detect the time evolution of a binary system. We 
work in natural units with G = c = h = I, but report all 
frequencies in Hertz. 



II. SPACE-BORNE DETECTORS 

A. Orbital effects 

The current design of the LISA mission calls for three 
identical spacecraft flying in an equilateral triangular for- 
mation about the Sun. The center of mass for the con- 
stellation, known as the guiding center, is in a circular 
orbit at I AU and 20° behind the Earth. In addition to 
the guiding center motion, the formation will cartwheel 
in a retrograde motion with a one year period (see Fig. 
^ . The detector motion introduces amplitude (AM) , fre- 
quency (FM) , and phase modulations (PM) into the grav- 
itational wave signals 0, 0| . The amplitude modula- 
tion is caused by the antenna pattern being swept across 
the sky. The phase modulation occurs when the differ- 
ing responses to the two gravitational wave polarizations 
are combined together. The frequency (Dopplcr) mod- 
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FIG. 1: The LISA mission configuration. The dashed line 
represents the orbit of the guiding center, which has a radius 
of f AU. 



ulation is due to the motion of the detector relative to 
the source. Since both the orbital and cartwheel motion 
have a period of one year, these modulations will show 
up as sidebands in the power spectrum separated from 
the instantaneous carrier frequency by integer values of 
the modulation frequency, = f/yr- 

To describe the coordinates of the detector we work in 
a heliocentric, ecliptic coordinate system. In this system 
the Sun is placed at the origin, the x-axis points in the 
direction of the vernal equinox, the x-axis is parallel to 
the orbital angular momentum vector of the Earth, and 
the y-axis is placed in the ecliptic to complete the right 
handed coordinate system. Ignoring the influence from 
other solar system bodies, the individual LISA spacecraft 
will follow independent Kcplerian orbits. The triangular 
formation comes about through the judicious selection of 
initial conditions. In Appendix A we derive the space- 
craft positions as a function of time. To second order in 
the eccentricity, the Cartesian coordinates of the space- 
craft are given by 

x{t) = i?cos(a) + ^ei?(^cos(2a-/3)-3cos(/3) 

+ -e2i? 3cos(3a-2/?)-f0 cos(q!) 
8 V 

-5cos(a-2/3)) 
y{t) = R sin(a) + ]^eR(^sin{2a - /?) - 3 sin(/3) 



f 



e^i? ( 3 sin(3a - 2/3) - 10 sin(a) 



-t-5 sin(a - 2/3) 
z{t) = — \/3ei? cos(a — /3) 

+VSe^R(^cos^{a- P) + 2sm^{a- (3) 



(1) 



In the above R — I AU is the radial distance to the 
guiding center, e is the eccentricity, a ~ 27r/,„t + k is the 
orbital phase of the guiding center, and /3 = 27rn/3 -I- A 
(n = 0, f , 2) is the relative phase of the spacecraft within 
the constellation. The parameters k and A give the initial 
ecliptic longitude and orientation of the constellation. 
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Using the above coordinates the instantaneous separa- 
tions between spacecraft are found to be 



15 sin I a — A H — 
6 



cos ^3(q; — A)^ 



Luit) = L[l + - 



15 sin (a — A ) 

V 6/ 



cos ^3(a — A)^ 



^23(0 = L(l 
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15 cos(a — A) 



+ cos (3{a — X)j 



(2) 



with L = 2y/3Re. From this it is seen that to Unear order 
in the eccentricity the detector arms are rigid. By setting 
the mean arm-length equal to those of the LISA baseline, 
L = 5 X 10^ m, the spacecraft orbits are found to have 
an eccentricity of e = 0.00965, which indicates that the 
second order effects are down by a factor of 100 relative 
to leading order. 



B. Gravitational Wave Description 

An arbitrary gravitational wave traveling in the k di- 
rection can be written as the linear sum of two indepen- 
dent polarization states. 



h(C) = /i+(e)e+ + /ix(e)e^ 



(3) 



where the wave variable ^ = i — fc • x gives the surfaces 
of constant phase. The polarization tensors are given by 



cos(27/')e^ — sm{2ip)e^ 
sin(2-!/')e+ -I- cos(2?/')e^ 



(4) 



where ip is the principle polarization angle and the basis 
tensors e"'' and are expressed in terms of two orthog- 
onal unit vectors, 



(5) 



These vectors, along with the propagation direction of 
the gravitational wave, form an orthonormal triad, which 
may be expressed as a function of the source location on 
the celestial sphere {0,(j)), 



u ~ cos 6 cos (j)x + cos 9 sin (j)y — 

V = sin (f>x — cos (p y 

k = — sin 9 cos <f>x — sin 9 sin (p y 



cos( 



(6) 



The above basis set is defined with respect to the 
barycenter reference frame. For a binary system - the 



standard gravitational wave source in the LISA band - it 
is natural to introduce another basis that is aligned with 
principle polarization axes, p and q, of the gravitational 
radiation. The orientation of the principle directions is 
chosen such that there is a 7r/2 phase delay between the 
two polarization states. The connection between the two 
basis sets is a rotation by the principle polarization angle 
Ip about the shared propagation direction k. 

We model the gravitational waves from a binary system 
according to 



(n) 



(7) 



where is the orbital phase. The instantaneous fre- 
quency of the n'^ gravitational wave harmonic is given 
by 



TT Ot 



(8) 



Unless the binary is highly eccentric or highly relativis- 
tic, the dominant emission will be quadrupolar, with fre- 
quency /(^) = /2(C) I and will be well described by the 
restricted post-Newtonian approximation: 



cos^ i) cos2*(^) 



fex(0-- ^ r, cos.sin2vI/(C). (9) 

Here M is the chirp mass, Dl is the luminosity distance 
and i is the inclination of the binary to the line of sight. 
Higher post-Newtonian corrections, eccentricity of the or- 
bit, and spin effects will introduce additional harmonics. 



III. DETECTOR RESPONSE: ANALYTICAL 

For two spatially separated test particles in free fall, 
the effect of a passing gravitational wave is to cause the 
proper distance between the masses to vary as a function 
of time. Finding the detector response reduces to solv- 
ing for the appropriate timelike and null geodesies in the 
spacetime with the metric 

ds^ = -{l+2(j))dt^ + {l-2(t)){dx'^+dy^+dz^)+h,jdx'dx^ . 

(10) 

In the above equation denotes the Newtonian poten- 
tial set up by various bodies in the Solar system and 
hij denotes the time-varying metric perturbation due 
to gravitational waves described in the previous section. 
The relevant geodesies arc those of the two spacecraft, 
xiiji), (T2), and the photons sent from spacecraft 1 to 
2, x^(A). We need to find the path taken by the photon 
that leaves spacecraft 1 at time ti and arrives at space- 
craft 2 at time t2, which amounts to a classic pursuit 
problem in curved spacetime. The calculation must take 
into account a host of factors, some due to the Newto- 
nian potential, and some due to the gravitational wave. 
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During the time taken for the photon to travel between 
the spacecraft, both effects are small and can be treated 
independently. 

The Newtonian potential leads to a variety of effects, 
such as a Shapiro time delay AL/L ~ Mq/R, gravi- 
tational redshift Avjv ^ MqL/E?^ deflection of light 
M ~ MqL/R^, and tidal flexing AL/L ~ MqL^/R^. 
Each of these effects is considerably larger than any of 
the effects caused by the passage of the gravitational 
wave, and they have to be subtracted before the grav- 
itational wave data analysis begins. The first step in 
the subtraction relies on us being able to accurately 
model the orbital phase shifts using the Solar System 
Ephemeris. The second step in the subtraction employs 
a high pass filter to remove the residuals from the or- 
bital fit, which occur at harmonics of the modulation fre- 
quency fm ~ 1/yi' — 3.2 X 10^^ Hz. The orbital effects, 
and the procedure for their removal, should be included 
in the full end-to-end model, even though they do not di- 
rectly affect the response of the detector to gravitational 
waves. 

The effect of the gravitational wave on the phase shift 
can be found by setting (/) = in Eq. H1U|) and solving the 
geodesic equation for the photons and the spacecraft in 
the metric perturbed by the gravitational wave. There 
are two equivalent approaches for finding the phase shift. 
The first approach is to find the Doppler shift of the pho- 
ton emitted by the first spacecraft and received by the 
second. The Doppler shift is then integrated with respect 
to time to give the phase shift. The Doppler derivation 
is given in Appendix B. The second approach is to inte- 
grate along the photons trajectory to find thepath length 
variation caused by the gravitational wave The ex- 
pressions given in Appendix B is valid to all orders in 
the spacecraft velocity w, and to first order in the gravita- 
tional wave strain h. However, as we explained in Ref. 0, 
it is hard to justify keeping terms of order vh given that 
V ^ 10~^. ft would take a phenomenally bright source, 
with a signal to noise ratio of ^ 10^, for the vh cross 
terms to be noticeable. Working to leading order in v 
and ft-, the path length variation for a photon propagat- 
ing from spacecraft i to spacecraft j is given by 

where fy(t) points from test mass i to mass j and h(^) is 
the gravitational wave tensor in the transversc-traceless 
gauge. The colon here denotes a double contraction, a : 
b = a'^h,j. 

Applying Eq. (|ll|l to a pair of orbiting spacecraft re- 
quires the careful evaluation of the fjj(f) unit vectors. 
This calculation is complicated by the motion of the 
spacecraft and the finite speed of light. For a photon 
emitted from spacecraft i at time ti and received at space- 
craft j at time tj the proper evaluation of the unit vectors 
is 



The distance the photon travels between spacecraft is 
given implicitly through the relationship 



(13) 



Here wc have used the fact that the reception time is 
the emission time plus the time of flight for the photon. 
We can numerically estimate the magnitude of this point 
ahead effect by expanding the photon propagation dis- 
tance in a u/c series: 

%(t,) = Ly(t.)(l + njiU) ■ v,(i,) + 0(i;2)) , (14) 

where "Vj{ti) is the velocity of spacecraft j and 

Lij{U) = ||xj(tj) - x,(ti)|| (15) 

is the instantaneous spacecraft separation. For the LISA 
mission with a mean arm-length of 5 x 10^ m and space- 
craft velocity v « 27r/,„i? sa 10~^, pointing ahead gives 
a first order effect of approximately 10^ m. For compari- 
son, the orbital effects given in Eq. (0) impart a variation 
in the photon propagation distance of 10^ m. 

An arbitrary gravitational wave can be decomposed 
into its frequency components: 



HI) 



df. 



(16) 



Tij {ti ) 



(12) 



Such a decomposition allows us to rewrite Eq. (|ll|l in 
the form 

/oo ^ 
D(/,t,fc) :h(/)e2-^«d/, (17) 
-OO 

where the one-arm detector tensor is given by 

D(/, t^k) = ^ (f,, (t) ® (O) ^(/, t, k) , (18) 
and the transfer function is 

TifXk) = sinc|^^(l-^%(t))^ 

xexp|^z^(l-fc.fy(t))^ . (19) 

Here /* = l/(27r^y ) is the transfer frequency for the 
ij-arm. The transfer functions arise from the interac- 
tion of the gravitational wave with the detector. For 
gravitational radiation whose frequency is greater than 
the transfer frequency the wave period is less than the 
light propagation time between spacecraft, which leads 
to a self-cancellation effect accounted for by the trans- 
fer functions. Below the transfer frequency the transfer 
functions approach unity. This leads to a natural divi- 
sion of the LISA bandwidth into high and low frequency 
regions, which will be exploited in a later section when 
we approximate the response of the detector. 
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The connection of Eq. Hll() to what is actually mea- 
sured depends on the design of the gravitational wave 
detector. The current proposal for LISA is to have 
each spacecraft measure two phases differences, one for 
each arm. The phase difference, as measured 

on spacecraft j, is found by comparing the phase of the 
received signal from spacecraft i against the outgoing sig- 
nal's phase that is traveling back to spacecraft i. Inherent 
in the phase difference measurements are both the grav- 
itational wave signal and noise contributions from laser 
phase noise C(t), shot noise n''{t), and acceleration noise 
n°(t): 

+n%iU) + S,,£it,)) . (20) 

Here the time ti is implicitly found through ti = tj — 
iij{ti). The subscripts on the noise components indicate 
the directional dependence of that component: Cy is the 
laser phase noise introduced by the laser on spacecraft j 
that is pointed toward spacecraft i, nfj is the shot noise 
in the photodetector on spacecraft j that is receiving a 
signal from spacecraft i, and is the projected accel- 
eration noise from the accelerometer on spacecraft j in 
the direction of spacecraft i. The position noise and path 
length variation arc converted into a phase difference by 
multiplying by the angular frequency of the laser, 27ri'o- 
Once the six phase differences are measured and 
telemetered down, the different interferometer signals 
can be synthesized. For example, the Michelson signal 
formed by using spacecraft 1 as the vertex craft is 

Si{t) = ^uihi) + $21 (i) - $13(^31) - $3iW , (21) 
where t2i and t^i are found from 

t31 = t-^3i(i3i). (22) 

However, due to the relatively large laser phase noise, the 
Michelson signal will not be a viable option. Instead a 
number of so called TDI signals will be used These 
signals are built by combining time-delayed Michelson 
signals in such a way as to reduce the overall laser phase 
noise down to a level that will not overwhelm the detec- 
tor's output. A particular example of a TDI variable is 
the X signal [l^: 

X{t) ^ $12(^21) + $2l(t) - $13(i3l) 

-$3l(t)-$12(i2l)-*2l(tl3) 

+ *13(t3l) + *3l(il2), (23) 

where the new times ii2, ii3, t'21, and t'^^i are defined 
through the implicit relationships 

ti2 = t-ii — ^12(^12) 
^13 ~ ^31 ~ ^13(^13) 
t'21 ~ ^13 — ^21(^21) 
^31 = ^12-41(4)- (24) 



By permutations of the indices similar forms for the Y 
and Z-signals can be constructed. 

By writing the response of the detector in a coordinate 
free manner we are able to apply this formalism to an ar- 
bitrary space-based mission. All that has to be changed 
are the spacecraft orbits. It should also be emphasized 
that the response is calculated entirely in the time do- 
main. In later sections we develop approximations to the 
full response by working in a hybrid time/frequency do- 
main. This hybrid approach assumes extra information 
about the sources, which allows us to develop explicit 
expressions for the detector response. 



IV. DETECTOR RESPONSE: NUMERICAL 
A. Noiseless response 

As an application of the equations presented in the 
previous section, we have simulated the response of the 
proposed LISA mission. The LISA Simulator [T3 | is de- 
signed to take an arbitrary gravitational waveform and 
output the full response of the detector. To apply the 
equations we have elected to work entirely in the helio- 
centric, ecliptic coordinate system. Therefore, all times 
are evaluated in terms of Solar System Barycentric (SSB) 
time. The conversion to the detector time is through the 
standard relationship dr — -y/l — v'^{t)dt, but since we 
only work to leading order in v the distinction is not 
made. (In practice there will be difficulties in synchro- 
nizing the clocks on the spacecraft |23|, but they do not 
trouble the simulations.) 

The positions of the spacecraft are calculated to sec- 
ond order in the eccentricity, Eq. which includes the 
leading order flexing motion of the array. Tidal effects, 
and third order terms in the eccentricity, arc neglected 
for now. 

One of the guaranteed sources for the LISA mission is 
the cataclysmic variable AM Canum Venaticorum. This 
binary star system is comprised of a low mass helium 
white dwarf that is transferring material to a more mas- 
sive white dwarf by way of Roche lobe overflow. AM 
CVn's orbital frequency of 0.972 mHz, and close prox- 
imity to the Earth (~100 pc) make it a good calibration 
binary for LISA. Shown in Fig. |21 is the simulated re- 
sponse to AM CVn expressed as a strain spectral density 
hf{f). Note that the barycenter gravitational wave sig- 
nal will be approximately monochromatic, however, the 
motion of LISA introduces modulations that cause the 
signal to spread over a range of frequencies . 

Another LISA source, but one whose event rate is 
poorly known, is the merger of two super-massive black 
holes. Figure 131 shows the simulated response of LISA to 
two 10^ Mq black holes coalescing at a redshift of z = 1. 
The observation tracks the final year before coalescence. 
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FIG. 2: The simulated X strain spectral density of AM CVn, 
demonstrating the induced modulations caused by the motion 
of the detector about the Sun. For reference, the dotted line 
is the average noise in this region of the spectrum. 
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FIG. 4: A realization of the Michelson noise for LISA, ex- 
pressed here as a strain spectral density. The dark line within 
the noise is a 128 bin rms value, while the rising curve is a 
standard LISA sensitivity curve. 
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FIG. 3: The Michelson response of LISA to two 10*' M© black 
holes coalescing aX z = 1 {Dl = 6.63 Gpc). The dotted line 
is the average Michelson noise in the detector. 



B. Noise 

Laser phase noise, photon shot noise, and acceleration 
noise are expected to be the dominant forms of noise 
in space-borne detectors. As previously discussed. Time 
Delay Interferometry is used to reduce the effects of the 
laser phase noise to a tolerable level. We assume that 
the TDI signal processing is properly implemented, and 
therefore neglect laser phase noise in our simulation. 

The simulation of the noise is done in the time do- 
main by drawing random numbers at each time step 
from a Gaussian distribution with unit variance and zero 
mean. For the white photon noise we then scale the ran- 
dom number by the shot noise spectral density defined 
in Ref. [H {Sp, = 1.0 x 10-^2 ni2/Hz). For the col- 



ored acceleration noise we begin by generating a white 
noise time series scaled by the acceleration noise spectral 
density (^acc = 9.0 x 10~^° m^/s^'/Hz), then integrate it 
twice to arrive at a colored time series. The integration 
introduces a falloff in the power spectrum that is 
characteristic of acceleration noise. The results of this 
procedure for the Michelson signal are shown in Fig. 0] 

Comparing this graph to a standard LISA sensitivity 
curve |24l |. a number of differences are apparent. The 
most obvious one is the lack of rise in the high frequency 
region. This is because the standard sensitivity curve 
folds the average detector response into the noise curve. 
The Sensitivity Curve Generator includes the all sky av- 
eraged and polarization averaged transfer function, which 
equals 3/5 at low frequencies and grows as above 
the transfer frequency. A secondary difference is in the 
overall normalization, as the Sensitivity Curve Cenerator 
scales the path length variations by the interferometer 
mean arm-length of L, while we scale the path length 
variations by the optical path length of 2L. 

To arrive at a simulation of the X noise we combine 
the noise elements as dictated by Eq. (|23|l . Doing so 
gives the results displayed in Fig.|Sl which agrees with the 
predicted results. To see this, we start with the analytical 
expression of the average Michelson noise curve shown in 
Fig.H 

1-1 <^*' 

which is derived in the appendix of Ref. |23 • In the above 
/* = 1/ (27rL) is the mean transfer frequency for an arm. 
Next, we note that the X signal is formed be differencing 
two Michelson signals, one time delayed by roughly twice 
the light travel time between spacecraft. Therefore, the 
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FIG. 5: A realization of the X noise for LISA, expressed here 
as a strain spectral density. The dark line is a 128 bin rms 
value for the noise. 



FIG. 6: A comparison of the simulated rms X noise (dotted) 
to the analytical prediction (solid). 



noise will enter in the X signal as 

nx{t) = nM{t) - UAiit - 2L) , 

which has a Fourier transform of 

nxif)^nMif){l-e-''f/f') 

and a power spectral density of 

.2f f 



Sxif)='isin'{j-]SM{f)- 



(26) 
(27) 

(28) 



The strain spectral density of the X noise is given by 



= 2 



sm — 



(/). 



(29) 



Shown in Fig. Elis a plot of hj (f) along with the aver- 
age from Fig. |S1 Although the derivation of the X noise 
strain spectral density assumed constant arm-lengths we 
see that there is excellent agreement between the pre- 
dicted results of Eq. H29|) and the simulation, which in- 
cluded the variations in the arms. 

Although Eqs. lUH), and ^ give the full re- 

sponse of a space-borne detector, they are analytically 
difficult to handle and time consuming to evaluate. For 
this reason we will now explore some approximations to 
the full response that use information about the input 
waveforms and a simplified description of the detector. 
These approximations not only aid in the development 
of data analysis techniques, but also give a greater in- 
sight into the workings of the detector. 



V. LOW FREQUENCY APPROXIMATION 



In sections Inl and Hill we saw that the full response of a 
space-borne gravitational wave detector was complicated 



by the intrinsic arm-length fluctuations, pointing ahead, 
and the signal-cancellation accounted for in the transfer 
functions. As a first approximation to the response of 
LISA we will neglect all of these effects. That is, we will 
work to linear order in the spacecraft positions, evaluate 
all spacecraft locations at a common time, and set the 
transfer functions to unity. It should be noted that this 
approximation was originally worked out by Cutler [l^ 
and can be viewed as an extension of the LIGO response 
to space-borne detectors. The transfer function T(/, t, k) 
can be set equal to unity when f ^ f^. For the LISA 
mission, whose bandwidth is 10~^ to 1 Hz, the transfer 
frequency has a mean value of /* = 0.00954 « 10~^ Hz. 

In the limit / ^ /* and /// ^ i the path length 
variation reduces to 



1 k ■ Tij (t) 

^ L[h,{t)®hAt)-Hm)] ■ 



(30) 



Working in terms of strains and neglecting noise, the 
Michelson signal from spacecraft 1 is given by 



6£i2{t ^ 2L) + Si2i{t - L) 
2L 

Siisit - 2L) + Shijt - L) 
2L 

2L 



(31) 



The last line follows from the condition f ^ f^. 

Using (jSJ, lU, and l|30|) the strain can be re-expressed 

as 



where 



Slit) = h+{^,{t))F+ + h^{^,{t))F> 

Ci{t) = t-k-^i{t) 

= t — R sin 9 cos (a{t) — cj)' 



(32) 



(33) 
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is the gravitational wave phase measured at spacecraft 1. 
The antenna beam pattern factors, F^{t) and F^{t), are 
given by 



1 



- [ sm{2iP)D+{t) + cos{2i')D'' {t) ) , (34) 



where 

D+it) 
D^t) 



Working to hnear order in the eccentricity, the Keplerian 
orbits given in yield 



ruit) ® fi2(t) ~ fisit) ® fi3(t)^ 

ri2(t) (g> rMt) - hsit) ® hsit)) : .(35) 



-36sin2(6') sin (2a(<) - 2A) 
(3 + cos(26')) ( cos(20) ^9 sin(2A) 



sin(4a(i) - 2A)j + sin(20) 
cos {4:a{t) - 2A) - 9cos(2A) 
-4\/3 sin(26') sin (3a(t) - 2A - 

-3sin(a(t)-2A + <; 



(36) 



and 



— \/3cos(6i)(9 cos (2A — 20) — cos (4a(t) 
-2A - 2(/))) - 6 sin(6') ( cos {Sa{t) - 2A - </> 



+3 cos (a(i)-2A + , 



(37) 



Equations l|32() to H37|l constitute the analytical formal- 
ism for the Low Frequency Approximation. These equa- 
tions are numerically quick to evaluate and can be han- 
dled analytically. As a point of reference, the strain pre- 
sented in Eq. (|32|l can be shown to be equivalent (most 
easily through a numerical comparison) to that derived 
by Cutler ^13]. 

To test the range of validity of this approximation we 
used The LISA Simulator (TLS) as a template to cal- 
culate the correlation between the full response and the 
Low Frequency Approximation (LFA), 



rif) 



\Stls|slfa; 



(38) 



TLS/ 



'LFA/ 



Using fixed random choices for the source location and 
orientation we systematically varied the gravitational 
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FIG. 7: The correlation between the Low Frequency Approx- 
imation and the full response of the LISA detector for a 
monochromatic source. The oscillatory structure at high fre- 
quencies is due to the transfer functions introduced by the 
full response. 



wave frequency and calculated the correlation at each 
frequency. The results of this calculation arc shown in 
Fig.CI 

We found that the Low Frequency Approximation has a 
strong correlation to the true response for frequencies be- 
low 3 mHz, at which point the correlation drops to 95%. 
The steep turn down in the correlation as the transfer 
frequency is approached is to be expected as the Low 
Frequency Approximation neglects the self-cancellation 
effects encoded in the transfer functions. The wiggles 
at higher frequencies are due to the transfer functions 
present in the full response template stls- The precise 
structure of these oscillations depends on the source lo- 
cation through the k ■ fij{t) dependence in the transfer 
functions. However, the turn down at 3 mHz is location 
independent. The location dependence does not become 
strongly evident until the correlation value has dropped 
to roughly zero. 

The significance of a particular correlation value is de- 
pendent on the signal-to-noise ratio of the source. For 
high S/N the effects neglected in the approximation will 
be detectable. Conversely, for a low S/N one may con- 
tinue to use the approximation at higher frequencies as 
the difference would not be noticeable. 



VI. RIGID ADIABATIC APPROXIMATION 

A. Response formalism 

The breakdown of the Low Frequency Approximation 
comes about through neglecting the transfer functions. 
As a second approximation to the LISA response we will 
now include the transfer functions, but continue to hold 
the detector rigid by working to leading order in the 
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spacecraft positions and evaluating all spacecraft loca- 
tions at the same instant of time. Such an approximation 
has been worked out before for the case of a stationary 
detector in , but here we extend it to include the 

motion of the detector. 

Physically this approximation can be viewed in the fol- 
lowing way. At an instant of time we hold the detector 
fixed and send photons up and back along the interfer- 
ometer arms and calculate the phase difference. We then 
increment the time by a small amount, moving the rigid 
detector to its new position in space, and repeat the pro- 
cess. This sequence of stationary states is the origin of 
the term "Adiabatic" for describing the approximation. 

For chirping sources the Adiabatic approximation re- 
quires that the frequency evolution / occurs on a 
timescale long compared to the light travel time in the 
interferometer: /// <C L. When this condition does 
not hold the Rigid Adiabatic Approximation is no longer 
valid and the full response should be used. In the limit 
f / f L the path length variation Hll|) reduces to 



T(/,t,fc) 



where the one-arm detector tensor is given by 

and the transfer function is 
,2/, 

/ 



(39) 



(40) 



1 



,it)) 



X exp 



2/. 



1 - fc • nj (t) 



(41) 



The Michelson signal is given by 

Sh2{t - 2L) + Shiit - L) 



Sl{t) 



2L 

Seisjt - 2L) + Sisijt - L) 
2L 



which may now be expressed as 

si(0 =5]D(fc,/„) :h„(e), 



(42) 



(43) 



and the round-trip detector tensor takes the form 

D(fc,/)-i 
and the round-trip transfer function is 



ia(g)a)T{a-kJ)-{b(g,b)T{b-kJ)] , (44) 



smc 



X exp 



-t-sinc 



X exp 



/(I -a-k) 

-i^(3 
2f^ 



a-k) 



fjl + a-k) 
2U 



k) 



(45) 
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FIG. 8; The correlation between the Rigid Adiabatic Approx- 
imation and the full response of the LISA detector. The turn 
down at ~500 mHz is due to neglecting the higher order ef- 
fects in the spacecraft positions. 



The time dependent unit vectors, d{t) and 6(t), are given 
by 



d{t) 

m 



X2(t)-Xi(t) 

L 

X3(t) -Xl(t) 

L 



(46) 



Collectively these equations are the analytical formalism 
for the Rigid Adiabatic Approximation. As with the Low 
Frequency Approximation^ the expressions are computa- 
tionally quick to evaluate and can be easily manipulated 
analytically. 

Figure |S1 shows the correlation between the full re- 
sponse and the Rigid Adiabatic Approximation for a 
monochromatic gravitational wave. Note that by includ- 
ing the transfer functions we are able to extend agree- 
ment with the full response two decades in frequency 
beyond where the Low Frequency Approximation broke 
down. The turn down at ^ 0.5 Hz comes about through 
neglecting the second order terms in the spacecraft po- 
sitions. As we described in Sec. Ill Al the second order 
orbital effects are down by two orders of magnitude in 
comparison to the linear order. This shows up in the 
Rigid Adiabatic Approximation through the transfer fre- 
quencies, which are evaluated for a rigid detector. Nor- 
mally the transfer frequencies are given by 



1 



(47) 



but for a rigid detector this reduces to the static form 
/* = \/{2ttL). The extension to higher orders in the 
orbital eccentricity can be done. The trade off is that 
the expressions become more complicated since the trans- 
fer frequencies would then become functions of time. In 
turn, this would require that each transfer frequency be 
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FIG. 9: A comparison of strain spectral densities as calculated 
using the Rigid Adiabatic Approximation (dashed) and the full 
response (solid) for a 5000Mq intermediate mass black hole 
binary at z = 1. For reference, the dotted line is the average 
Michelson noise in the detector. 



FIG. 10: The correlation between a monochromatic Rigid 
Adiabatic Approximation and one that includes a 2PN phase 
evolution. The solid line represents a WD- WD binary with 
mass components of 0.5 Mq, the dashed line is for a NS-NS 
binary with masses 1.4 Mq, and the dotted line is for a 10 
Mq black hole and a 1.4 M© neutron star binary. 



evaluated along each arm during each time step rather 
than using one constant value throughout the entire cal- 
culation. Additionally, the normalization of the unit 
vectors in Eq. would need to be evaluated at each 
step since the arm-lengths would vary as a function of 
time via Eq. (01 • Such an approach would be appropri- 
ately called the Flexing Adiabatic Approximation since 
the arm-lengths would now oscillate in time about a mean 
value of L. Although the expressions would become an- 
alytically complicated, the numerical evaluation would 
not be significantly slower since the additional steps are 
straightforward to evaluate. 

Figure ini compares the output of The LISA Simulator 
to the Rigid Adiabatic Approximation for a binary system 
of intermediate mass black holes, each of mass 5000Mq 
at a redshift of z = 1. The observation covers the final 
year before coalescence. The agreement is excellent. 

In a recent paper Scto used a variant of the 

Rigid Adiabatic Approximation to calculate the effects of 
LISA'S finite arm-lengths on the analysis of gravitational 
waves from chirping supermassive binary black holes. A 
comparison of our Rigid Adiabatic Approximation, which 
is derived from the path length variation, and Seto's ap- 
proach, which is based on the Doppler shift formalism is 
given in Appendix C. 



B. Applications 

Utilizing the speed of the Rigid Adiabatic Approxima- 
tion we may investigate various data analysis questions. 
Here we provide one concrete example by determining 
when phase evolution of a binary system due to radiation 
reaction needs to be included in the source modeling. 

For our calculations we used the restricted Post- 



Newtonian approximation, whereby the gravitational 
wave amplitude is calculated to first order, while the 
phase evolution is calculated to second order The 
justification for this is that LISA will be far more sen- 
sitive to the phase than the amphtude The lack 
of additional harmonics of the orbital period also simpli- 
fies the calculation as we only have to calculate a single 
transfer function at each time step. 

To quantify the importance of including the evolution 
of the gravitational wave phase, we calculated the corre- 
lation between a monochromatic Rigid Adiabatic Approx- 
imation to one in which the phase evolution is included. 
Figure ^1 shows the correlation for three types of bina- 
ries expected to reside inside our own galaxy: a white 
dwarf binary with mass components 0.5 M0, a neutron 
star binary with masses 1.4 M0, and a 10 Mq black hole 
with a neutron star companion. 

What we found is that the frequency at which the 
monochromatic signal diverges from one that includes 
phase evolution depends on the masses of the binary com- 
ponents. The reason for this comes from the expression 
for /, which contains a mass dependent coefficient. For 
the stellar mass binaries we studied, the drop in the cor- 
relation happened to coincide with the breakdown of the 
Low Frequency Approximation. Thus, most Milky Way 
sources can be modeled as monochromatic sources using 
the Low Frequency Approximation. 

Another way to represent the same data is to map the 
initial frequency to the time of coalescence. The results of 
this calculation are shown in Fig. In this case we see 
that chirping becomes important for stellar mass sources 
within ~10^ years of coalescence. As expected, the map- 
ping to the new variable preserves the mass dependence 
seen in Fig. ^1 
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FIG. 11: The correlation between a monochromatic Rigid 
Adiabatic Approximation and one that includes a 2PN phase 
evolution, this time expressed in terms of the time to coales- 
cence. 

A final way to represent this data is to set the indepen- 
dent variable equal to the change in the frequency scaled 
by a bin width, 6f — {ff — fi)//S.J, where for one year 
of observation the bin width is A/ = 1/yr ~ 3.2 x 10~^ 
Hz. Such an approach is shown in Fig.^l Unlike with 
the other representations of the correlation between a 
monochromatic and coalescing signal, the results of this 
calculation are independent of the system's masses. It is 
also interesting to note that this result implies that it will 
be possible to detect if a source is coalescing or not well 
within a bin width. This fact is not in conflict with the 
Nyquist theorem, which states that the frequency resolu- 
tion will not be better than the inverse of the observation 
time. The reason being that we have additional informa- 
tion, namely the functional form of the phase evolution, 
which is not assumed in deriving Nyquist 's theorem. 
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FIG. 12: The correlation as a function of the fractional bin 
width change in the frequency, 5f = (// — fi)/Af, between a 
monochromatic Rigid Adiabatic Approximation and one that 
includes a 2PN phase evolution. Note that all three types of 
sources previously considered are included in this plot. 

information matrices in studies of astrophysical parame- 
ter extraction. On the other hand. The LISA Simulator 
is available if we need to simulate the response to highly 
relativistic gravitational wave sources such as the merger 
of two black holes. 
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APPENDIX A: SPACECRAFT POSITIONS 



VII. DISCUSSION 

We have examined the forward modelling of space- 
borne gravitational wave detectors with special emphasis 
on the LISA observatory. Forward modelling will play 
two distinct roles in the developement of space-borne ob- 
servatories. The first is as part of a complete end-to-end 
model that takes into account every concievable physi- 
cal effect, and the second is as an intermdeiary between 
source simulation and data analysis. Here we have fo- 
cusscd on the latter role, and to that end we have stud- 
ied two simple approximations to the full response - the 
Low Frequency Approximation and the Rigid Adiabatic 
Approximation. We found that the Rigid Adiabatic Ap- 
proximation could be used in place of the full response 
for a wide range of data analysis projects. For example, 
the relatively simple analytic form of the Rigid Adiabatic 
Approximation is well suited to the calculation of Fisher 



For a constellation of spacecraft in individual Keple- 
rian orbits with an inclination of i = -\/3e the coordinates 
of each spacecraft are given by the expressions 

X = r |^cos(-\/3e) cos f3 cos 7 — sin /3 sin 7^ 

y = r ^cos(A/3e) sin/3cos7 -f- cos/3sin7^ 

z = — r sin(v^e) C0S7 (Al) 

where (3 = 2n7r/3 -I- A (n = 0, 1, 2) is the relative orbital 
phase of each spacecraft in the constellation, 7 is the 
ecliptic longitude, and r is the standard Keplerian radius 



1 + e cos 7 



(A2) 



Here R is the semi-major axis of the guiding center and 
has an approximate value of 1 AU. 
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FIG. 13: The angle k, gives the initial ecliptic longitude of the 
guiding center. 




positions of the spacecraft as functions of time, 
x{t) = R cos(q;) + ii?e ^ cos(2q! — /3) — 3 cos(/3) 

+ lRe^ (3 cos(3a - 2/3) - 10 cos(a) 



-5cos(a-2/3) 

y{t) = Rsm{a) + ii?e(sin(2a - (3) - 3sin(/3) 

+ lRe^ (3 sin(3a - 2/3) - 10 sin(a) 
8 V 

+5sin(a-2/3)) 

z(t) = — \/3i?ecos(a — /?) 

+\/3i?e^(cos2(a-/3) + 2sin2(a-/3)^ .(A8) 

These are the desired coordinates of each spacecraft as a 
function of time. 



APPENDIX B: DOPPLER STYLE DERIVATION 
OF THE FULL RESPONSE 



FIG. 14: As viewed by an observer at the origin, A gives the 
initial orientation of the spacecraft constellation. 



To get the above coordinates as a function of time we 
first note that the ecliptic longitude is related to the ec- 
centric anomaly, ip, by 



tan 



1 + e /V 

tan — 

1 - e V 2 



(A3) 



and the eccentric anomaly is related to the orbital phase 
a{t) = 2iTt/T + K via Kepler's equation 



(A4) 



a ~ 13 ^ ip ~ e sin ip . 



Assuming a small eccentricity we may solve Eq. (IA4|) 
through an iterative process where we treat the e sin tp 
term as being lower order than -0, 



(A5) 



= a - /3 + e sin ipn-i ■ 

Through such a procedure we arrive at 

^ — a — (3 + e sin{a — l3) + e^ cos(q; — /?) sin(Q; — /?) + ■•■ . 

,(A6) 

Substituting this result into Eq. (jA3p and expanding to 
second order in the eccentricity gives an ecliptic longitude 
of 



7 = (a-/3)-|-2esin(a-/3) + -e^ cos(a-/3) sin(a-/3)H . 

(A7) 

Substituting the ecliptic longitude series into Eq. (|Aip 
and keeping terms up to order gives the Cartesian 



The Doppler shift of a photon emitted by spacecraft 1 
and received by spacecraft 2 can be elegantly derived 
using the symmetries of the spacetime (|10|l . When (/) = 
the spacetime admits three Killing vectors, 



C(i) = C( 



2) 



C(3) = 5t + fc . 



(Bl) 



These provide three constants of the motion, (^(jj • [/, 
which along with the normalization condition [/•[/ = 
or J7- J7 = — 1, fully specify U{\) in terms of some initial 
four velocity U (0). Writing the metric as g^^ = 77^1^ 

y 1 

we may express the photon propagation four-vector as 
= - l^rh^ps'' , (B2) 

where s is a null vector in the unperturbed geometry: 
s^s'^rjf^i, — 0. At the time of emission from spacecraft 1, 
s{ti) = So + Ssi, while at the time of reception at space- 
craft 2, s(t2) = So + 6s2. Here so is parallel to the unit 
vector connecting the two spacecraft in the unperturbed 
spacetime, while Ssi and Ss2 are perturbations to the 
path due to lensing by the gravitational wave. Defining 



As" = s"(i2)-s"(ii). 



we have 



1 



fT"(t2) = a"(ti) + As" - -r^'^i'Ahp^s 

which yields 

2s"As'3r;„^ = cj^{t2)cj^{t2)9c.i3(t2) 

= 0. 



(B3) 



(B4) 



(B5) 
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Equations ljBl|l . ljB3|) and (|B5p yield four equations for 
the four As": 

s'^AsPr^^p - . (B6) 
These can be solved to give, for example, 



2k ■ . 



(B7) 



Here k — > (l,fc) is the null propagation vector for the 
gravitational wave. The frequencies of the emitted and 
received photons, as measured at spacecraft 1 and space- 
craft 2, are given by vi = ~Ui{ti) ■ (f{ti) and 1/2 = 
—U2{t2) ■(^{^2) respectively. Here Ui and U2 are the four- 
velocities of the two spacecraft. Note the vi = is the 
operating frequency of the laser on board spacecraft 1. 
Evaluating v ~ —U-a yields 

V = -7(t) (^s\t) + v\{t)s^{t)f^,, + hi,^{ty{t)sl 

(B8) 

where = j = dt/dr and the are the ordinary three 
velocities of the spacecraft. The spacecraft trajectories 
U may be expressed in terms of the unperturbed trajec- 
tories Uo according to 



(B9) 



where A" are constants set by the initial conditions at 
some time t. Once the initial conditions for the space- 
craft have been set, equations (jB6p . (jB8|) and (jB9p give 
the full Doppler shift 1^2 — vi at any subsequent time. 
The expressions simplify considerably if we drop terms 
of order , vh and higher: 



Aiy = V2 — vi ~ As* 



(BIO) 



2k ■ So 

Converting this into a fractional frequency shift, Av/vq, 
and using Sq = i/q a, where d is the unit vector connecting 
the two spacecraft in the background geometry, we have 

Aiy d®d: Ah /t^hn 

(Bll) 



vo 2(1 - fc • a) 

Integrating the above expression with respect to time 
yields the time delay described by Eq. ((TT|l . 



APPENDIX C: RECONCILIATION BETWEEN 
ALTERNATIVE RIGID ADIBATIC 
FORMALISMS 

According to Eq. (2.2) of Ref. |2i|, the relative fre- 
quency shift for a photon traveling from spacecraft 2 to 
spacecraft 1 can be expressed as 

ysiit) = ^(^+cos(2V'i2) + iAxSin(2i/'i2)) 

X (1 - cos 012) {U{t, 1) - U{t - L, 2)) ,(C1) 



where the function U{t,i) gives the phase of the grav- 
itational wave at spacecraft z, 0ij is the angle between 
the source location on the sky, —k, and the detector arm 

X.; - X , 



cos ( 



k ■ fji 



and tpij is given through the relationship 



tan-0i 



rj^ ■ p 



(C2) 



(C3) 



Here p and q are unit vectors along the principle polar- 
ization axes of the gravitational wave. The amplitude 
coefficients, and Ax, are given as functions of the 
orbital inclination angle and the intrinsic amplitude A, 



A+ = ^(1 + cos^ l) 
Ax ~ 2A cos i . 



(C4) 



By direct substitution of Eqs. (|(]2p . H('3p . and l|('4p into 
the relative frequency shift we arrive at 



ysiit) = 



{hi ■ pf + {f2i ■ qf 



x(A+((f21-p)2-(f21-g)2) 

+2iAx{f2i ■p){f 21 •?)) 
(C/(i,l)-C/(t-L,2)) 



2(l-fc-f2l) 



The overall coefBcicnt 



l-(A:-f2i)2 



(r2i • pf + {hi ■ qf ' 



(C5) 



(C6) 



can be shown to equal unity by writing the vector in 
terms of the orthonormal triad {p, q, fc} according to 



fij = sin 6' cos (t)'p + sin 6' sin 0'^ + cos 9 k . 



(C7) 



To proceed further we note that the gravitational wave 
basis tensors can be expressed as functions of p and 



p ® p — q® q 
p® q + q® p . 



(C8) 



Using these relationships, it follows that 

{hj ® hj) ■■ e+ = {hj ■ pf - (fy • qf (C9) 
(fy (g) fij) : = 2{hj ■ p){fij ■ q) . (CIO) 



Returning to the expression for the relative frequency 
shift we now have 



2/31 (0 



'"21 (Xl ^21 



2(l-fc-f2i) 
x(C/(i,l)- [/(t-L,2)) 
^21 ® f21 ■ Ah 
2(l-fc-f2i) ' 



iAxe' 



(Cll) 
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where in the last step we combined the amphtudc and wave evaluated at each spacecraft. The above results 
phase functions to form the difference in the gravitational agrees with Eq. (|ljll() found in Appendix B. 
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